#### cronbach's alpha

cronbaa1 <- function(id_vec = rep(1:5, each = 4), 
                     answer_vec = sample(1:4, length(id_vec), replace = T),
                     ques_vec = rep(1:4, 5)
                     ){
  k <- length(unique(ques_vec))
  
  sums <- c()
  for(i in unique(id_vec)){
    sums <- c(sums, sum(answer_vec[id_vec == i]))
  }
  tvar <- var(sums)
  
  add_var <- 0
  for(i in unique(ques_vec)){
    var_temp <- var(answer_vec[ques_vec == i])
    add_var <- add_var + var_temp
  }
  alpha <- (k/(k - 1)) * (1 - (add_var/tvar))
  return(alpha)
}


## ---
# res <- c()
# for(i in 1:100){
#   res <- c(res, cronbaa1())
# }
# hist(res, breaks = 40)

# id_vec <- rep(1:50, each = 4)
# answer_vec <- rep(c(4,4,3,4,
#                 1,2,3,1,
#                 1,1,2,1,
#                 3,3,4,4,
#                 2,3,3,1),10)
# ques_vec = rep(1:4, 5)
# 
# cronbaa1(id_vec, answer_vec, ques_vec)
